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We present a new approach to the problem of binary black holes in the pre-coalescence stage, 
i.e. when the notion of orbit has still some meaning. Contrary to previous numerical treatments 
which are based on the initial value formulation of general relativity on a (3-dimensional) spacelike 
hypersurface, our approach deals with the full (4-dimensional) spacetime. This permits a rigorous 
definition of the orbital angular velocity. Neglecting the gravitational radiation reaction, we assume 
that the black holes move on closed circular orbits, which amounts to endowing the spacetime with 
a helical Killing vector. We discuss the choice of the spacetime manifold, the desired properties 
of the spacetime metric, as well as the choice of the rotation state for the black holes. As a 
I simplifying assumption, the space 3-metric is approximated by a conformally flat one. The problem 

is then reduced to solving five of the ten Einstein equations, which are derived here, as well as the 
boundary conditions on the black hole surfaces and at spatial infinity. We exhibit the remaining five 
^ , Einstein equations and propose to use them to evaluate the error induced by the conformal flatness 

_ O ■ approximation. The orbital angular velocity of the system is computed through a requirement which 

' reduces to the classical virial theorem at the Newtonian limit. 
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I. BACKGROUND AND MOTIVATION 



~^ • Binary black holes have been the subject of numerous studies in the past two decades, both from the analytical and 
^ numerical point of view. These studies are motivated by the fact that the coalescence of two black holes is expected 
O to be one of the strongest sources of gravitational waves detectable by the interferometric detectors LIGO, GEO600, 

TAMA300 and VIRGO, currently coming on-line ||. 
^ From the analytical point of view, the most recent works are based on the post-Newtonian formalism (see e.g. 

Ref. [|j for a review) or on the effective one-body approach developed by Buonanno and Damour |^|-^. In these 
^ works, the black holes are treated as point mass particles^], which is a very good approximation when the black holes 
are far apart. For closer configurations, one may turn instead to some numerical approach. The numerical studies 
can be divided in two classes: (i) the initial value problem for two black holes (see Ref. for a review) and (ii) 
the time evolution of the initial data (see Ref. [|| for a review and Refs. for recent results). One of the major 

problems in this respect is to get physically relevant initial data. Indeed, initial data representing two black holes have 
been obtained long ago by Misner and Lindquist jl^], as well as Brill and Lindquist |jl^ (a modern discussion of 
these solutions can be found in Ref' [|l5[ or Appendices A and B of Ref. [|l6|). However these solutions correspond to 
two momentarily static black holes and are therefore far from representing some stage in the evolution of an isolated 
binary black hole in our universe. 

Based on the seminal work of Bowen and York |1^, Kulkarni, Shepley and York |Q have described a procedure to 
get initial data representing binary black hole with arbitrary positions, masses, spins and momenta. This procedure, 
known as conformal imaging, has been used by Cook and other authors ||ig|-[2^ to numerically construct 3-dimensional 
spacelike hypersurfaces representing these initial data. These solutions constitute some generalization to the non-static 
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regime of the Misner-Lindquist solution ||T^,|l3| , the spacehke hypersurface having the same topology: two isometric 
asymptotically-flat sheets connected by two Einstein- Rosen bridges. Also based on the Bowen and York's work 
another approach, the so-called puncture method, has been undertaken by Brandt and Briigmann and used 
recently by Baumgarte The resulting solutions also have arbitrary positions, masses, spins and momenta but 

their topology is different from that of the conformal-imaging approach: the spacelike hypersurface has now three 
asymptotically-flat sheets, which are not isometric. These three sheets are connected among themselves by two 
Einstein-Rosen bridges; one of the sheets contains two throats and is supposed to represent our universe. These 
solutions hence constitute some generalization to the non-static case of the Brill-Lindquist solution |0. Both the 
conformal-imaging and the puncture approaches assume that the metric of the spacelike hypersurface is conformally 
flat. A third approach to the problem relaxes this assumption; it has been developed recently by Matzner, Maroronetti 
and collaborators p^-|2^. These authors use a linear combination of two boosted Kerr-Schild metrics as the conformal 
3-metric in York's treatment of initial conditions . 

The main drawback of the three approaches described above is that they contain freely specifiable parameters, 
related to the values of positions, momenta and spins of the black holes, and it is difficult to figure out which 
parameters correspond to a physical configuration. In particular, it is not obvious at all how to select, among all these 
configurations, those that correspond to binary black holes on closed circular orbits. Of course, the circular orbits are 
approximate representations of the exact orbital motion, which is inspiralling due to the loss of energy and angular 
momentum via gravitational radiation. However in the regime where the radiation reaction time scale is much longer 
than the orbital time scale, i.e. before the last stable orbit, we expect that the binary system can be approximated 
by a sequence of tighter and tighter circular orbits. Note that the gravitational radiation reaction makes initially 
elliptic orbits become circula r [pO[ . It is thus legitimate to search for such orbits. This problem has been addressed by 
Cook [ pT[[_ and Pfeiffer et al. ||22| within the conformal-imaging approach and by Baumgarte within the puncture 
approjichg (see e.g. for a review of these computations). Although they differ in the topology of the spacelike 
hypersurface (two-sheeted for |2^j2^ against three-sheeted for both sets of studies rely on the effective-potential 

method proposed by Cook . This method amounts to defining the binding energy by a somewhat ad- hoc formula, 
and to define the angular velocity of a circular orbit by minimizing the binding energy with respect to the angular 
momentum at fixed total energy and separation. This method can be criticized on the following ground: the only well 
defined global quantities on an asymptotically flat spacelike hypersurface are the ADM total mass M and ADM total 
linear momentum. The latter can be chosen to be zero without any loss of generality. With some restrictions on the 
asymptotic gauge, the total angular momentum J can also be defined ||3^ . Defining the binding energy would require 
the notion of individual mass for each hole, let say Mi and M2, in order to set E'bind = M — Mi — M2. However there 
does not exist any unique definition of the individual masses Mi and M2 for a binary black hole. In Refs. |pT| , p2| , p5| 
the authors define Mi and M2 via the formula which relates the mass of a Kerr black hole to its horizon area and 
its angular momentum. However such a formula is strictly demonstrated only for a isolated rotating black hole (Kerr 
spacetime). In particular, it does not take into account any tidal effect. 

Our approach for finding circular orbits of binary black holes is very different: instead of considering 3-dimensional 
spacelike hypersurfaces, we adopt from the very beginning a 4-dimensional point of view, i.e. we consider a full 
spacetime containing two moving black holes. Of course, in order to make the problem tractable, we introduce some 
approximations, the most significant being the assumption of strictly circular orbits, which amounts to endowing our 
spacetime with a Killing vector field (helical symmetry). An interesting pay-off is that the orbital angular velocity 

can be unambiguously defined as the rotation rate of the Killing field with respect to some asymptotically inertial 
observers. This definition does not suffer from the ambiguity of the 3-dimensional approaches and is made possible 
only because we have re-introduced time in the problem. 

The presentation of this new approach is organized as follows. We set up the problem in Sec. ||, starting from the 
explicit construction of the spacetime manifold, introducing a metric, as well as a corresponding isometry on it, and 
finally imposing the helical symmetry. The Einstein equations are then considered in Sec. p^ , first in their general 
form and then after the assumption of a conformally flat 3-metric. Global quantities such as the ADM mass and 
the total angular momentum are discussed in this section, as well as the virial prescription for the orbital velocity. 



Section tV deals with the asymptotic behavior of the shift vector and the extrinsic curvature tensor and discusses the 
connection between helical symmetry and asymptotic flatness. Finally Sec. contains the concluding remarks. 



■^To our knowledge, the Kerr-Schild approach has not been used yet to get circular orbits, the article providing results 
only for black holes in hyperbolic motions. 
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II. FORMULATION 



A. Spacetime manifold 



1. Construction 



We consider the s pac etime to be a difFerentiable manifold with the topology of the real line M times the Misner- 
Lindquist manifold More precisely, for any couple of positive numbers (ai, 02) and any couple of real numbers 

{xi,X2) such that |a;i — 2;2| > oi + 02, let us consider the subset of M'^ obtained by removing the interior of balls of 
radius ai and 02 and center x = xi and x = X2'- 

£ := {{x,y,z) €W^,{x - xif +y'^ + >ala.nA {x - X2f + >al] . (1) 

Let us call Si and 5*2 the 2-spheres defining the "inner" boundaries of £: 

Si:^{{x,y,z)£W',[x-xif+y'' + z'^ai] , (2) 



S2:^{ix,y,z)emMx-X2f+y^ + z^ = al} . (3) 

Let us consider two copies £1 and £11 of £ and define Aii :=M.x £1 and Aiu := R x £u. The spacetime manifold Ai is 
be then defined as the union A4i U A^n with both and S2 of each copy identified (see Fig. |^). The reader is referred 
to Sect. IV of Ref. or Sect. II of Ref. for a precise construction of the manifold structure in the vicinity of Si 
and ^2. The part A^i of AI will be designed hereafter as the upper space and the part A^n as the lower space. The 
boundaries iSi := M x 5*1 and 1S2 := K x <S'2 between Mi and Aiu are called respectively throat 1 and throat 2. 




FIG. 1. Construction of the spacetime manifold. 
Hereafter we label by (i, xi, yi, 01) the points of M.i considered as a part of M 



M? {£1 being a part of R^), and 
The corresponding two charts will be called 
the canonical coordinate systems. These two charts cover Ai minus the two throats. The whole manifold Ai can be 
covered entirely by a single coordinate system: 



by (t, xii, yii, zii) the points of Aiii considered as a part of 



Ci : Ai 
P 



{t,xi,yi,zi) if P e All 
Ii{t, xii, yii, zii) if P e AAii 



(4) 



where Ii 



Ji(t,x,y,z 



^ denotes the inversion through the 2-sphcre Si 
a\(x — X\) 



a\y 



[x — xiY' + + ' (a; — xi)^ +7/^ + 2^' (x ~ x\)^ + y'^ + z^ 

In a similar way, one can introduce the coordinate system C2 associated with throat 2. 



(5) 
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In the coordinate system Ci or C2, the throats are not located at constant coordinate values. Therefore, it is more 
convenient to introduce instead the polar coordinate system (t,ri,6i,ipi) centered on throat 1, as follows: 



for P e Mi, 



and 



for P e Mil, 



xi — ri sin 9i cos (pi + xi 
Ui — ri sin 61 sin (fi 
zi = ri cos 9i 



a? 

xii — ^ sin 9i cos ipi 
yii = ^ sin 6*1 sin ipi 
zii = 77 COS 01 



- Xi 



(ai < ri < +00) 



(6) 



(0 < ri < ai) 



(7) 



The throat 1 corresponds to ri = ai. The polar coordinate system {t,r2, 02,^2) centered on throat 2 is introduced 
similarly. Note that any of the two coordinate systems {t,ri,6i,Lpi) and {t,r2, 62,^2) covers the whole spacetime 
manifold M. For the {t,ri,9i,ipi) system, Mi corresponds to ai < ri < +00 and Mii to < ri < ai. Similarly, for 
(t, r2, 02,^2) system, Mi corresponds to 02 < ^2 < +00 and Mii to < r2 < 02 (see Fig. ||). 























( 

























FIG. 2. Coordinate systems {t, ri,9i,ipi) and {t, r^, O2, ^2) on the spacetime manifold M. Shown here is a t = const section 
of M, with the dimension in the 6 direction suppressed, leaving only (ri, ipi) or (r2, 952). 



2. Canonical mapping 

From the very construction of M , we have at our disposal the canonical mapping (see Fig. |^) 

/: Mi — > Mil 

{t,xi,yi, zi) I — > {t = t,xii = xi,yii ^ yi, zii ^ zi) 



(8) 



Note that in terms of the (t, ri,9i,ipi) coordinate system, this map can be written as an inversion through the sphere 
ri = ai (see Fig. |): 



l{t,ri,0i,ipi) = [t,^,0i,ipi 



In terms of the {t, r2, 02,^2) coordinate system, it looks like an inversion as well: 



I{t,r2,e2,ip2) = t, ^,02,^2 

r2 



(9) 



(10) 
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Let (x") be a coordinate system on A4i [for instance {x") = {t,xi,yi, zi), (x") = (t,ri,9i,ipi) or {x") = 
{t,r2, 62,^2)] and (j/^) be a coordinate system on Mu [for instance (y'^) = (t, xn, yn, 2:11), (y'') = {t,ri,6i,ipi) 
or (y'') = (i, r2, 92,^2)]- The map / is fully characterized by its components with respect to the coordinates (x") 
and (y''): the image I{P) of a point P £ A^i with coordinates {x'^ , x^ , x'^ , x^) has the coordinates 

y^ = /''(x°,a;\x2,a;3) . (11) 

Let us now examine the action of the map / on various fields on A^. If / is a scalar field on Mu, I induces a scalar 
field on Ml through 

hf:=foI, (12) 

i.e. 

hf[p] = mp)] , (13) 

for any point P of Mi- 

Also / maps any vector field v on TWi to a vector field /*v on M\\ through 

/,v(/) :=v(/,/) , (14) 

for any scalar field / on Mn- If vectors are represented by their components with respect to the coordinate bases 
djdx'^ and d/dy^, one has 



and, according to definitions ( |14| ) and (|12D, 

/,v(/) = — /(/'^(x'^)) . (16) 

Hence the matrix of the mapping between vectors on A4i and vectors on A4u is given by the Jacobian matrix of /: 

(W[/(P)] = |^z;"[P], (17) 

where P denotes any point of A^i. 

The action of / on vectors can be used to define the action of / on 1-forms as follows: / maps any 1-form a; on Mu 
to a 1-form I^u) on Mi through 

:= uj{hM) , (18) 

for any vector field v on Mi- If 1-forms are represented by their components with respect to the coordinate bases 
dx°' and dyf^, one has 

IM^) = (/,a;)ada;"(v) = (/*a;)„7;" (19) 

and, according to definition (p^), 

where the third equality arises from Eq. (^7|). Comparing Eqs. (|l^) and (|2^) leads to 

(/,a;)4P] = |^a;,[/(P)] (21) 

at any point P in Mi. 

Similarly, the action of / on bilinear forms can be defined as follows: / associates any bilinear form T on Mu to a 
bilinear form on Mi according to 

/*T(v,w) :=T(/,v,/,w) . (22) 
One can show easily that in terms of the components with respect to the coordinate bases dx" (S> dx^ and dy^ (g) dy'', 

(^*T).,[^] = |^|5^-['(^)] (23) 

at any point P in Mi- 
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B. Spacetime metric 



1. Properties 

We endow Ai with a Lorentzian metric g with the following properties: 
(1) g is asymptotically fiat at the ends of Mi and Mu : 



, g-^, (24) 

, lim g = »7 , (25) 
^11+^11+^11^°° 



where rj is a flat metric. 

(2) the canonical mapping / is an isometry of g: 

/*g = g • (26) 

(3) The t — const sections of A4 are maximal spacelike hypersurfaces with respect to g. 

The assumption (1) is introduced because w e cons ider only isolated systems. Its connection with the quasi- 



stationarity hypothesis will be discussed in Sees. [IC 1 and IV C. 

The assumption (2) is motivated by the fact that the Schwarzschild and Kerr spacetimes possess such an isom- 
etry. This can be readily seen when using isotropic (quasi-isotropic for Kerr) coordinates instead of the standard 
Schwarzschild (Boyer-Lindquist) ones. By virtue of Eq. ([23|), the isometry condition (^) can be expressed in terms 
of the components of g at any point P in A^i : 

|^5^5.^W^)]=5«M^]. (27) 

It is also useful to write the isometry condition on the contravariant components of the metric tensor; by means of a 
generalization of Eq. (p7|), one gets: 

The assumption (3) is motivated by the well-known good properties of maximal slicing among which there 

is the singularity avoidance. 



2. 3+1 decomposition 

In this article we use the standard 3-1-1 formalism of general relativity foliating the spacetime by a family of 
spacelike hypersurfaces. From the very construction of A4, a, natural foliation is by the t = const hypersurfaces E(, 
where t is the same coordinate as that introduced above. By virtue of assumption (3), this constitutes a maximal 
slicing of spacetime. 

Let us denote by n the future directed unit normal to St. Being normal to St, n should be coUinear to the gradient 
of i: 

n = -NVt , (29) 

where N is the lapse function, which can be seen as a normalization factor such to ensure that n • n := g(n, n) = —1. 
Let us now examine the behavior of n under the isometry I. By the definition (^, / preserves the hypersurface St. 
According to the definition (^^, the square of the norm of /,n is 

g(/,n,/,n) = (/,g)(n,n) . (30) 

But thanks to the isometry condition /*g = g, the last term in this equation is simply n • n = — 1. Hence 

/*n • /*n = -1 , (31) 

i.e. / preserves the norm of n. Similarly, for any vector v tangent to St, / preserves the scalar product n • v = 0, so 
that /*n • J,v = 0. But since Sf is globally invariant under /, /*v represents any vector tangent to St, so that /*n is 
in fact normal to Sf. Having the same norm than n, we conclude that 
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(32) 



Since t, considered as a scalar field on A^, is preserved by /, so is its gradient and the relation 
( ^ ) results then in the following transformation law for the lapse function: 

LN = ±N . 



9h, combined with 



(33) 



In order to understand the significance of the ± sign in Eqs. (^2|) and (p3[), let us consider the case of a single 
static black hole, i.e. the (extended) Schwarzschild spacetime. Two kinds of maximal slicing of this spacetime are 
depicted in a Kruskal diagram in Fig. ^ starting from the same initial hypersurface v ~ 0, t ~ 0. The first one 
corresponds to a symmetric lapse [sign + in Eqs. ( |33| ) and (^2[)]. The throat is located at it = 0; the slicing penetrates 
under the event horizon {R = 2M), and accumulates on the spacelike hypersurface R = 1.5M |^,^. The second 
slicing corresponds to an antisymmetric lapse [sign — in Eqs. ( p^ ) and (p2|)]. In fact, it corresponds to the standard 
Schwarzschild solution in isotropic coordinates: 



M 

1 + — ) [dr^ + r^{de^ + sin^ Odip^)] 



(34) 



with 



N : 



1 - M/2r 
1 + M/2r 



(35) 



This lapse function is clearly antisymmetric under the transformation I : r ^ K'P /{Ar) across the throat located at 
r — M/2. The negative value of the lapse for r < M/2 is easily understandable when looking to Fig. while t is 
running upward on the right part of the diagram (corresponding to A^i), it is running downward on the left part 
(corresponding to Adu). Since in the Kruskal diagram the future direction is everywhere upward, the lapse should be 
negative in A4u. 




FIG. 3. Kruskal diagrams showing the slicing of the extended Schwarzschild spacetime by two families of maximal hyper- 
surfaces; left: lapse function symmetric with respect to the isometry I; right: lapse function antisymmetric with respect to the 
isometry I. R denotes the standard Schwarzschild radial coordinate and r the isotropic one. 

Let us now consider a coordinate system (x') on each Sj. For instance, it can be chosen in one of the three 
coordinate atlas introduced so far: {(xi, yi, zj), {xu,yu, zu)}, {(ri, 6*1, (/3i)} and {(^2, 6*2, </32)}- (i,a;') constitutes then 
a coordinate system on Ai. The shift vector (3 associated with the coordinates (i,a;') is defined by the following 
orthogonal split of the coordinate basis vector d/dt: 

d 

— = iVn + /3 with n • ^ = . (36) 
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Since the transformation / is purely spatial d/dt is preserved by it. By virtue of Eqs. ( p2| ) and (JS^), the product Nn 
is also invariant with respect to /. Consequently 

= /3 . (37) 

The S-metric induced by g on the hypersurfaces is 

7 g + n (8) n . (38) 

From Eqs. ( ^ ) and (^2|), we obtain immediately that / is also an isometry for the 3-metric 7: 

/*7 - 7 ■ (39) 

The components of the metric tensor can be expressed in terms of the lapse function and the components of the 
shift vector and the 3-metric, according to the classical formula 

gf^^dx^'dx" = -{N^ - f3^P') dt^ + 2/3, dt dx' + 7,^- dx^ dx^ . (40) 

The extrinsic curvature tensor K of the hypersurface St is given by the Lie derivative of the 3-metric along the flow 
defined by the normal to Eii 

K = -i£n7- (41) 
By the symmetry properties ( |3^ ) and (^9|), we obtain that 

/,K = ±K. (42) 



3. Explicit isometry conditions in polar coordinates 

In what follows, we consider o nly po lar coordinate systems centered on one of the two throats, i.e. either the system 
{t, ri,6i,(fi) introduced in Sec. 11 A 1 or {t, r2, 02, <f2)- For the sake of clarity we will drop the indices 1 or 2 on r, 6 
and if. It should be understood that the formulas will be valid for cither coordinate system. The Jacobian matrix of 
/ with respect to {t, r, 9, cp) is easily deduced from Eq. (^ or ([lO|): 

5^=diag(^l, 1, ij , (43) 

where a denotes either a\ or 02. From the isometry condition ( p8| ) expressed on g", we get, since 17" — —l/N"^ : 

N[I{P)f = N[P]^ , (44) 

for any point point P in A4i, i.e. we recover the already established relation (|3^. From the isometry condition ( [2^ ) 
expressed on g** we get 

/3'-[I{P)] = -^P^[P] (45) 

P'[I{P)]^(3%P] (46) 
/3^[/(P)] =/3^[P] , (47) 

where we have used g'' — (3^ /N'^ and Eq. (^) to go from g** to /3'. Again note that we recover the isometry condition 

®.- . .. ^ 

Finally the isometry condition (E^) expressed on ~ results in 

lrr[P] = '^Irrim] (48) 
"1r0[P\ = ~^lr0[I{P)] (49) 
lr^[P\ = -^lr^[I{P)] (50) 

l0e[P]=le9[I{P)] (51) 
le^[P] = IB^P)] (52) 
l^AP] = l^Am] ■ (53) 
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Comparing Eqs. ( p6| ) and (^2|), we see that the isometries properties of the components Kij of K are the same as 
those above for 7.^ , except possibly for an opposite sign. 



4- Choice of the isometry sign 

As discussed above, the behavior of the foHation with respect to the isometry / involves a + or — sign in the 
transformation rules of the unit normal [Eq. (|3^)], lapse function [Eq. (|33|)] and extrinsic curvature [Eq. (j^)]. In 
this article, we choose the sign to be the minus one. This is motivated by the fact that the maximal slicing with the 
+ sign of the Schwarzschild spacetime (left part of Fig. |3|) does not respect the stationarity of the problem, i.e. the 
Killing vector djdt of Schwarzschild geometry does not carry a slice of that foliation into another slice |Q (see also 
Sec. IV of Ref. (S^]). On the contrary, the slicing with the — sign (right part of Fig. ^ respects the stationarity of the 
problem. It seems to us more appealing to use a slicing which in the case of a single black hole, makes the problem 
time- independent. We regard the artificial time dependence resulting from the + sign as an unnecessary complication. 
Beside simplicity, another advantage of the — sign choice is to allow us to test the numerical code by comparison with 
the standard form of the Schwarzschild or Kerr metric in the special case of a single black hole. 

Thus, from now on, we set 



(54) 



LN = -N 



and 



hK = K . 

Eq. (^) can be explicited for any point P in A4i: 

N[I{P)] = ~N[P] , 
which amounts to choosing the — sign when taking the square root of Eq. 



(55) 
(56) 
(57) 



I- 



5. Boundary conditions on the throats 

An immediate consequence of Eq. (^7|) is that the lapse function vanishes on the two throats: 

N\s, = and N\s^ = . 



(58) 



Indeed from the very definition of / [Eq. (g)] and the construction of A4 by identifications of the two copies of Si or 
^2, every point P in Si or 52 is a fixed point for /. Hence Eq. ( |57| ) results in N[P] = ~N[P] on Si and S2- 
Similarly, Eq. (^) implies that the r component of the shift vector vanishes on the throats: 



and 



(59) 



Taking the first derivatives of Eqs. ( 
of the isometry of the shift vector: 



)-(47), we get the additional following relations on the throats, as a consequence 



39 

dip 

dr 
8(3^ 
dr 



= 
= 
= 
= 



(60) 
(61) 
(62) 
(63) 
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where we have dropped the indices 1 or 2 on r, 0, v? and S. Note that relations (|6^) and (^l]) could have been obtained 
also as consequences of Eq. (^9|) since the throats are located at a constant value of the coordinate r. 
Equations (E8h-(M) and their first derivatives give the following values for the 3-nietric on the throats: 



dr 










djre 




— 


80 


s 

















de 


s 










dr 


S 




dr 


s 










dr 


s 









and 

and 



. 



dip 

d^r^p 



dip 



(64) 
(65) 
(66) 
(67) 
(68) 

(69) 

(70) 

(71) 



6. Apparent horizons 

As a direct consequence of the isometry hypothesis, the throats Si and ^2 are minimal 2-surfaces of the spatial 
hypersurface Sf. Moreover, as shown by Cook & York the fact that K is antisymmetric with respect to the 
isometry / [Eq. (|56|)] implies that Si and ^2 are apparent horizons. 



C. Quasi-stationarity hypothesis 



1. Helical Killing vector 



As discussed in Sect. |, we consider binary black holes in the quasi-steady stage, i.e. prior to any orbital instability, 
so that the notion of closed circular orbits is meaningful. Following Detweiler [Q, we translate these assumptions 
in terms of the spacetime geometry by demanding that there exists a Killing vector field £ such that, near spacelike 
infinity. 



_d_ 

dfo 



d_ 

dipa 



(72) 



where to and ipa are respectively the time coordinate and the azimuthal coordinate associated with an asymptotically 
inertial observer, and O is a constant, representing the orbital angular velocity with respect to the asymptotically 
inertial observer. Let us call £ the helical Killing vector. We refer the reader to for a detailed description of this 
concept. 

The helical symmetry amounts to neglecting outgoing gravitational radiation in the dynamics of spacetime. For 
non-axisymmetric systems — as binaries are — it is well k nown that imposing £ as an exact Killing vector leads 
to a spacetime which is not asymptotically flat In Sec. IV C , we will exhibit explicitly how the deviation from 
asymptotical flatness arises. However, from a physical point of view, the exact helical symmetry is too strong an 
assumption because it assumes that the binary system is rotating on a fixed orbit since the past time infinity. Doing 
so, it has filled the entire space with gravitational waves, such that their total energy is a diverging quantity, whence the 
impossibility of asymptotic flatness. A weaker assumption, which is compatible with asymptotic flatness and sounds 
physically more reasonable, is the following one. Due to the reaction to gravitational radiation the binary system is 
in fact spiraling. Therefore in the past time infinity, it was infinitely separated. As a consequence, the amount of 
emitted gravitational waves was very weak. The integral of their energy density is now a converging quantity, which 



10 



allows for asymptotic flatness. The quasi-stationarity hypothesis should then be understood as imposing a helical 
Killing vector on a part of spacetime limited in time. 

It is natural to demand that the isometry associated with the Killing vector £ preserves, not only {Ai, g) as a whole, 
but also the sub-structure of A4 defined by M-i, A4ii and the two throats. This amounts to demanding that for any 
of the coordinates system {t, a:') introduced above, where t is the coordinate used explicitly in the construction of Ai, 

The above equality means that t is an ignorable coordinate. It does not mean that the problem is stationary in the 
usual sense of this word, for t is not timelike at spatial infinity: by virtue of relation ([72[), i ■ i ^ ^(x\ + y\) > 
when xi, yi — > 00. 

2. Rotation states of ttie black holes 

The above geometrical assumptions are intended to correspond to a physical system of two black holes in a quasi- 
steady state. We have not specified yet the rotation state (spin) of each black hole. In this article, we consider 
synchronized (or corotating) black holes. This rotation state can be translated geometrically by demanding that 



the two throats be Killing horizons |41| associated with the helical symmetry. This means that each null-geodesic 
generator of Si and 52 must be parallel to £. In particular, this implies that the Killing vector £ is a null vector on 
the throats: 

£ ■ £\g^ =0 and £■ e.\s^ = . (74) 

As a guideline, note that this condition is verified by the helical Killing vector d/dto + flud/dtpo of the Kerr space- 
time, where d/dto, d/dipo and fin are respectively the Killing vector associated with stationarity, the Killing vector 
associated with axisymmetry and the rotation angular velocity of the black hole. This classical result is known as the 
rigidity theorem in the black hole literature ]4^ ] . 

In a recent work, Friedman et al. |^9| note that corotation (in the above sense) is the only possible rotation state 
consistent with the helical symmetry in the full Einstein theory. However, a weaker definition of quasi-equilibrium 
(not assuming that £ is an exact Killing vector, as we do here) allows for more general rotation states, as shown very 
recently by Cook 143]. 

Combining Eqs. ([zs]) and ( |36| ) shows that £. is related to the lapse function, unit normal and shift vector through 

£ = Nn + l3 , (75) 

so that the scalar square of £ is 

£-£ = -N'^+l3-l3 . (76) 

Thanks to the vanishing of the lapse on the throats, the rigidity condition ( |7^ ) is then equivalent to /3 ■ /3 = on iSi 
and ^2. But /3 being a vector parallel to St, /3 • /3 = 7(/3,/3); the positive definitcncss of the 3-metric 7 implies then 

f3\s^=0 and I3\s,^0. (77) 

Hence, not only the r-component of (3 is zero [Eq. (|59|)], but the total vector f3 vanishes on the throats. 

III. EINSTEIN EQUATIONS 
A. General form 

The vacuum Einstein equations can be written as the Hamiltonian constraint equation: 

R - K,jK'^ = , (78) 

the momentum constraint equation: 

DjK'^ = , (79) 
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and the "dynamical" equations: 

dK 

-^~£f^K,, = -D,DjN + N{R,,-2K,kK)) , (80) 

where Rij denotes the Ricci tensor of the S-metric j, R = R^ the Ricci curvature scalar, and Di the covariant 
derivative associated with 7. Note th at we have used the vanishing of the trace of K, as a consequence of the maximal 
slicing (assumption (3) in Sec. II Bl). Besides, the geometrical relation ( pl| ) involving the extrinsic curvature results 
in the following equation: 

^ - £(31-^3 = -'^NK,, . (81) 
Following York , Shibata & Nakamura , and Baumgarte & Shapiro , we introduce the "conformal metric" 

:= -i-^'H, , (82) 

where 7 is the determinant of the 3-metric components 7ij. 7^^ is a tensor density of weight —2/3. York |Q] has 
shown that it carries the dynamics of the gravitational field. One can introduce on Ef a covariant derivative Di such 
that 

(i) A7y = ; 

(ii) if 7ij is conformally flat (7^ — "if^fij), then Di = Di, where Di is the covariant derivative associated with the 
flat metric fij. 

We refer to Refs. for details in the case of Cartesian coordinates and to Ref. for any coordinate system. 

Note that the property (i) is not sufficient to fully characterize Di since the covariant derivative Di fulfills it as well, 
reflecting the fact that 7^^ is a metric density and not a metric tensor: there exists at least two distinct covariant 
derivatives "associated" with it. Let us denote by Rij the Ricci tensor associated with the covariant derivative Di 
and by R the corresponding scalar density: R := '^^''Rku where 7'-' is the inverse conformal metric 

:= 7i/3^''^' . (83) 

Let us also introduce the following tensor densities 

A,j := j-^/^K,j and A'^ := j^/^K'^ , (84) 

and denote by D^ the operator ^^'^Df^. The Hamiltonian constraint equation (^) can then be written as an equation 
for the determinant 7: 

A^Mn7+^Aln7^'ln7= ^ (^-7'/'^^'') • (85) 
The momentum constraint equation ( |79| ) becomes 

DjA'^ + ^A'^Dj ln7 = . (86) 
The dynamical Einstein equations ( pO[ ) can be decomposed into their trace part 

DiD'N + ln7L>'iV = -fiNAijA'^ (87) 

and their traceless part 

N (f''j''Rki + —& In 7 A ln7 ) + - (d' ln7 AiV + A ln7 Aiv) - -y^^/^DW^ {-f^^'^N) 
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iv(^i?+^Aln7^''ln7) +^Aln7^'-'^- 7~'^^^fc^''(7'/'^) 
+7I/3 (^2N%iA'''A^' + p'^DkA'^ - A^Wkp' - A^'^DuP^ + '^DuP^A'^^ = 
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Note that in Eqs. (|8^) and (|8^), we have used the helical symmetry to set to zero the time derivatives and that we 
have explicited the Lie derivatives along (3. Similarly, the evolution equation (|8l| ) for jij can be split into its trace 
part 



and its traceless part 



A/3^ = -^/3'Aln7 (89) 



2NA'^ = + b^P' - \bkl3^ . (90) 



Inserting this relation into the momentum constraint (g^) results in the following equation for the shift vector: 

DjD^f]' + ^DWjP^ + f^Rjkp'' + A'^ (^NDj In 7 - 2DjN^ = . (91) 

We recognize here the minimal distortion equation of Smarr & York [^3|, i.e. we recover the fact that the shift vector 
of coordinates co-moving with respect to a Killing vector field is necessarily a minimal distortion shift. 

B. Approximation of a conformally flat 3-metric 

1. Equations 

As a first step in this research project, we introduce the approximation of a conformally flat 3-metric: 

7 = , (92) 
^ being some scalar field, a nd f the canonical flat 3-metric associated with the canonical coordinates {xi,yi,zi) and 



(xii, yii, zii) (see Sec. P A 1| ) 



Such an approximation has been used in all previous studies of binary black hole initial data based on the conformal 
imaging approach |1{^-P3|1 or on the puncture approach p4J23]. It has been relaxed in the recently developed Kerr- 



Schild approach |2£~|2^. Strictly speaking, the assumption (|9^ ) is exact only for a single non-rotating (Schwarzschild) 
black hole. However, as discussed by Mathews et al. [Q, such an approximation is quite good even for a maximally 
rotating Kerr black hole. 

As an immediate consequence of Eq. (b3) , we have 



7 = *'V, (93) 
where / is the determinant of the metric components fij . The conformal "metric" takes then the simple form 

l^J=r'^'f^J and f^=f/^fi. (94) 



By property (ii) of the covariant derivative Di (see Sect. Ill A above), another consequence of Eq. ( |92| ) is that 

b, = A , (95) 
where bi is the covariant derivative associated with the flat metric f. Note that, from the definition of D*, one has 

b' = /i/^L)' , (96) 

where l)* := p^bk- It follows immediately from Eq. ( p5| ) that the Ricci tensor Rij is identically zero. 
The relations ( p^ can be rewritten 

A,, = r^'^A,, and A'^ = f^'^A'^ , (97) 

where we have introduced the tensor fields: 

iy := -^^^K.j and A'^ := 'f'^K'^ . (98) 
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Taking into account the above relations, the Hamiltonian constraint equation ([Sq) becomes an elliptic equation for 



A* 



where A := D}^D^ is the Laplacian operator with respect to the flat metric f. 
The momentum constraint equation, under the form becomes 

A/3* + \&Dj(3^ 2i*^ (DjN - 6NDj In*) , 

whereas the trace part of the dynamical Einstein equations, Eq. (|87|), becomes 

AiV = N^'^AijA'i - 2D J In ^ D^ N . 
The traceless dynamical Einstein equations (|8|) reduces to 



(99) 



(100) 



(101) 



2ND' In * D^ In * + In D' N + D^ In * IJ'iV - -^''^ D' DH-^"^ N) 

4 



1 

3 

vl/4 



2NDk In ^Z?'' In * + 2Dk In ^'D'^iV - -*"2A(*2iV) 



2NfkiA'''Ai^ + /3''Dk{A'^) - A'^Wkp' - A'^Dkfi' + -Dkl3^ A'^ 



= 



The trace part of the evolution equation for , Eq. ( p9| ) , becomes 

A/3* = -6/3* A In* , 

whereas its traceless part ( pO| ) results in a relation between the extrinsic curvature tensor and the shift vector: 
where (LPY^ denotes the (flat) conformal Killing operator Es] applied to the vector /3: 



(102) 



(103) 



(104) 



(L/3)* 



D'pi + A/3* - -Dkf3'' 
3 



(105) 



2. Solution scheme 



Our approach is the following one: using Eq. (104) to evaluate A'^ , consider Eqs. (|9S|), ( |100| ) and ([101^ as coupled 
elliptic equations to be solved for respectively *, (3 and N. The remaining five Einstein equations, Eqs. ( |102| ), are not 
used to get the solution. Moreover, they are not satisfied by the solution (*, (3, N), except in special circumstances 
(e.g. spherical symmetry). This reflects the fact that the conformally flat form ( |92| ) constitutes only an approximation 
to the exact Einstein equations. An interesting application of Eqs. ( |102| ) is then to evaluate its left-hand side in order 
to gauge the error resulting from the conformally flat approximation. Besides, note that Eq. (103) is not used in the 
above scheme. We will disc uss t his point in Sec. IIIB4. 



The system of Eqs. (|98|)-( |l0l[ ), resulting from the assumption of conformal flatness and maximal slicing, has been 
already proposed by Isenberg & Nester , as well as Wilson & Mathews ||5l[ , as an interesting approximation to the 
full Einstein equations. It has been used by many authors to compute binary neutron stars on circular orbits [p2|-p9[. 



3. Boundary conditions 



The equations (|99|), (100) and (101) we are facing being elliptic, it is very important to discuss the boundary 
conditions to set on their solutions. Thanks to the isometry /, the computational domain is chosen to be half the 
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full spacetime, i.e. only A4i. Its boundaries are then the spatial infinit y and the two throats Si and S2- At 
infinity, the metric should be asymptotically flat (hypothesis (1) of Sec. II B 1 ). This implies that 



1 when ri ^ oo or r2 



and 



— > 1 when ri — > oo or 7-2 — > oo . 
Combining Eqs. ([72|), ( fz^ ) and (107), we get the asymptotic behavior of the shift vector: 



/3 



17—^ when ri 

dtpo 



00 or r2 



The boundary conditions on the throats have been derived in Sec. II B, In particular, Eqs. (|6J), 
equivalent to the the following condition on the conformal factor 'i>: 







\dri 





= 



and 



dr2 2r2 



= . 



52 



and 



All the remaining equations listed in Eqs. (|64|)-([71|) are automatically satisfied by the conformally flat form 
boundary equation on the lapse have already been given [Eq. (p8|)]: 

N\s, = and N\s^ = . 

as well as that on the shift vector, resulting from the rigid rotation hypothesis [Eq. (|77|)]: 

(3\s^ = and (3\s^ - . 



spatial 
(106) 

(107) 

(108) 
1]) are 

(109) 
, The 

(110) 
(111) 



4-. Regularity on the throats and isometry of the shift 



A direct consequence of Eqs. ( |llO| ) and (104) is that the shift vector on the throats should satisfy not only (111) 
but also 







and 



152 



, 



(112) 



in order for the extrinsic curvature to be regular on the throat. Note that in the case of a single rotating black hole, 
such a condition is equivalent to 8(3'^ /dr = and dp^ jdQ — 0. The first condition follows from e.g. Eq. (10.25) of 
Ref. |4|] and the last one from the rigidity theorem (Jj"^ is constant - and zero - on the horizon). In the present case, 
the properties ( p9| ) and (|60|)-(|6^), which follow from the isometry /, in conjunction with the property (111), which 
follows from the rigidity assumption R, imply that 



and 



(L/3)-|< 



- -1 diy 



3 dr 

2 9/3'' 

2 9/3'- 



(L/3) 



re I 



ir'^sva^e dr 



s 



(113) 

(114) 

(115) 

(116) 
(117) 



^Note that the constancy of /3 on the throat, implied by Eq. (Ill), results in the vanishing of all the partial derivatives of the 
components /3' with respect to and if. 
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Therefore the condition (112) is equivalent to 



dri 



= 



and 



5i 



= . 



(118) 



Now the trace of th e re lation between the extrinsic curvature and the derivative of the 3-metric, Eq. (103), gives, 
when inserting Eq. (Ill) in its right-hand side. 



l5i 







and 







(119) 



From Eq. (|113| ), it follows that Eq. Qll^ ) is satisfi ed a s well. This establishes the regularity property ( |112D . 

The above argument relies on the fact that Eq. (103), relating the trace of the extrinsic curvature tensor (here zero ) 
to the divergence of the shift vector, is satisfied. However, as discussed in Sec. [II B 2| , we solve only Eqs. (|9^)-( |l0i| ) 
to get the metric fiel ds ^ , N and (3. This means that there is a priori no guarantee that Eq. (103) is satisfied by 
the solution of (^-( |lOll ) (see Sec. IV. C of Ref. |39| for a discussion of this point). It has been argued recently by 
Cook that if one reformulates the problem by assuming that the helical vector £ is not an exact Killing vector, 
but only an appro xima te one — as it is in reality — then the on ly fr eely specifiable part of the extrinsic curvature, 
as initial data, is (104), not (103). This means that the relation (|l03| ) between the extrinsic curvature and the shift 



is not as robust as the relation (104). 

Another problem is that, when solving the system (|9S|)-(101) subject to the boundary conditions (|l09| )-( lll ), there 
is no guarantee not only that the solution for the shift vector obeys to Eg. (|103|) but also that it obeys to the isometry 
conditions (|6|) and ( 63|) [the other isometry conditions, namely (^9|), (|60[) a!nd (|6ll), are satisfied by virtue of the 
boun dary condition (111)]. In the companion article | |60[ |, we p resent a method to enforce the regularity condition 
(118) as well as the isometry condition conditions (|6|) and (|63|). This amounts to add, at each step of the iteration, 
a corrective term f3^^^ to the solution /3constr of the momentum constraint (100), so that the shift vector 



/3 = /3c 



/3c 



(120) 



is well behaved, i.e. satisfies (i) the rigidity boundary conditions (111), (ii) the condition (11?) which ensures the 
regularity of the extrinsic curvature on the throats, and (iii) the isometry condition (|3^). 

If at the end of the iteration, /3(.or converged to zero, then we get a regular solution of the Einstein equations 
in the conformal flatness approximation. On the contrary, if P^^^. stays at some finite value, we get a solution which 
violates the momentum constraint equation. The numerical solutions we have computed p0[ | belong to this category. 
However they have (cf. Sec. IV. B of Q) 



\f3,J < 10-3|/3| 



(121) 



which shows that the momentum constraint is only very slightly violated. Taking into account the other approxima- 
tions performed, e.g. conformal flatness, we find this to be very satisfactory. 



C. Global quantities 

The total mass-energy content in a St hypersurface is given by the Arnowitt-Deser-Misner (ADM) mass M, which 
is expressed by means of the surface integral at spatial infinity 

M = ^ r V {Dj^ki - Dki.i) dS, (122) 

(see e.g. Eq. (20.9) of Ref. ^). In the case of the conformally flat 3-metric 7.y — '^^fij, this integral can be written 

M = — - I D'^^dS^ . (123) 
27r Joe 

By means of the Green-Ostrogradski formula, this expression can be converted into the volume integral of A^* plus 
surface int egra ls on the throats; using the Hamiltonian constraint ( p9|) to express A^P, as well as the boundary 
condition (109) on the throats, one gets 



M = — I xfjd^x + — / * sin6li dOi dwi + — I ^ sm02 d92 d(p2 ■ 

167r7 "^T^ Jr^=a^ '^^Jv,=a2 



(124) 
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The total angular momentum in a Ej hypersurface is defined by the following surface integral at spatial infinity 



J = 



T? dSi 



(125) 



where m 

asymptotic). Note that in the present case K\ = (maximal slicing). Note also that J defined by Eq. (125) coincides 



d/difiQ [see Eqs. ( [72[ ) and (108)] is the rotational Killing vector of the flat metric f (to which 7 is 



with the z-component of t he v ector J' defined by Eq. (12) of Bo wcn & York |17[]. As discussed by York ||29| , |32| , 
contrary to the definition ( |122D of the ADM mass, the definition (125) of J requires some asymptotic gauge-fixing 
conditions stronger than the mere asymptotic flatness (p^)-(p5|), because of the supertranslations ambiguity. Some 
natural gauge-fixing conditions are provided by the asymptotic quasi-isotropic gauge proposed by York |29[| . Such 
conditions are fulfilled by the conformally flat metric (^). Using the fact that \1/ = 1 at spatial infinity, we can 
replace K'^^m^ by ^'^A^-'/jfem'^ in Eq. ( |125| ) and express J by means of the Green- Ostrogradski formula as a volume 
integral plus surface integrals on the throats. The volume integral vanishes identically, as one can see by considering 
the following identity: 



A(*'i'V,fem'=) = A(*'i'^)/,fem''^ + ivI;6i«J[5,(/,,m'=) +5,(/,fem'=)] . 



The first term on the r.h.s. vanishes by virtue of the momentum constraint (p6D, which can be written as 



(126) 



(127) 



whereas the second term vanishes for m is a Killing vector of the flat metric f. Thus the formula for J reduces to 
integrals on the throats: 



OTT Jri=ai ^TT Jr2=a2 



(128) 



where dSi denotes the surface element with respect to the flat metric f and oriented toward the "interior" of the 
throats. 



D. Determination of the orbital velocity 



The orbital angular velocity does not appear in the partial differential equations listed in Sec. IIIB 1. It shows 



up only in the boundary condition (108) for the shift vector. This contrasts with the binary neutron star case, where 
enters in the equation governing the equilibrium of the fluid (see e.g. ^T\). 

At this point, it appears that, solving Eqs. (|99|)-(101), with the boundary conditions ( |l06[ )-( p^8| ) and (109)-(111), 
one c an get a solution {N,f3,'^) for any given value of fl. For instance, if we set = in the boundary condition 
( |108D , we get /3 = as a solution of (IOC) and the Misner-Lindquist solution for |lj,|l^. Of course, such a solution 
is not admissible on physical grounds, and we need an extra condition to fix fi. 

As a boundary condition at spat ial i n finity , we have demanded only that g tends to the Minkowski metric of flat 
spacetime [conditions (24)-(p5|) or ( |l07| )-( p^ )]. We could go a little further and demand instead that g tends to the 
Schwarzschild metric corresponding to the ADM mass M. This implies the following behavior for the conformal factor 
^ and the lapse N [cf. Eqs. dH) and (H)]: 



^ - 1 



— when r 
2r 



(129) 



N - 1 



M 
r 



when 



(130) 



where r denotes either the coordinate ri or r2. From the very definition of M, the behavior (12£) is guaranteed by 
Eq. (|l2|). However, the solution of ( p^ ) is such that iV - 1 - A/'/r, with a priori M' ^ M. The behavior ( p^ ) is 
thus a extra condition imposed on the system (p9|)-(101). This is the condition which will enable us to fix VL. 

Let us show that for stationary spacetimes, i.e. i n the case where t is timelike at infinity, the extra condition (|130| ) 
follows from the remaining Einstein equations (102), i.e. the equations that we have not used in the system (99)- ([101^. 
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Indeed, the quadratic terms of the type In ^D^N or In In ^ which appear in Eq. ( 102 ) all decay at least as 



wh en r — > oo. Now, for stationary spaceti mes, it can be seen that the L ie derivative along (3 of A^^ which appear 



in Eq. ( |102|) , decays also at least as r^** [Eq. ( |197|) below]. Then Eq. (|102|) implies that DW^^-^^N) decays at least 



as r , which means that the l/r (monopolar) part of N vanishes, i.e. 



— when r 



(131) 



This is possible only if and N have opposite monopolar l/r terms, which implies the property ( |130| ). 

Note that, for stationary spacetimes , the monopolar term of the lapse N is the Komar mass associated with the 
timelike Killing vector. The condition ( |l30| ) is then intimately linked to the virial theorem: as already shown by two 
of us , a relativistic generalization of t he cl assical virial theorem can be obtained provided that the Komar mass 
coincides with the ADM mass [property ( |130[ )]. This last property has been shown to hold for asymptotically flat 
stationary spa cetim es by Beig [ |6^ . In order to exhibit more clearly the link with the virial theorem, let us combine 
Eqs. (pi) and (|lOl]) to derive an equation for ^^A^ (see e.g. Eq. (51) of Ref. E^): 



2Di^D\'^N) 



where we have re-introduced a non-vanishing stress-energy tensor T^^, via Saf3 



discussion when considering the Newtonian limit. The condition (131) is equivalent to the vanishing of the monopolar 
term of ^I^^A^, i.e. from Eq. (132) and assuming a spacelike slice Si diffeomorphic to R'^, 



(132) 



Tau for the benefit of the 



AnS, 



-Ai,A'^ 
4 ^ 



2 A* D'i^N) \ ^/f(fx = 



(133) 



It is easy to see that this relation is equivalent to the relativistic virial theorem giv en b y Eq. (29) of Ref. [£2j , once 
the latter is specialized to a conformally flat 3- metric. The Newtonian limit of Eq. (133) is nothing but the classical 
virial theorem: 



2T + 3P + W = 



(134) 



where T is the total kinetic energy of the system, P the volume integral of the pressure and W the gravitational 
potential energy. Note that th e val ue of for two Newtonian particles of individual mass m in circular orbit (radius 
R) can be obtained from Eq. (134) {T = mR^Vl^, P = Q, W ^ — m^/(2i?)); this results of course in the Keplerian 
value ^2 = 2m/{2Rf. 

Let us remark that Detweiler [Q has proposed to determine the orbital velocity of binary black holes in circular 
orbits by means of a var iatio nal principle. Although he does not state it precisely, his variational principle also use 
the "virial" assumption (130) [cf. the not so well justified sentence "In the gauge described in Chapter 19 of Misner 
et al. (1973) the flux integral at infinity is 47rM - SttJ" below his Eq. (17)]. 



E. Generalized Smarr formula 



A formula relating the ADM mass M , the total angular momentum J, the angular velocity O and some integrals 
on the throats can be obtained as follows. The key point is to notice that the Einstein equations ( [ts]), ([79|), the trace 
of ( |80| ) and (^ij) imply, when d/dt is a Killing vector [Eq. ([73|)], the following remarkable identity [|j8[ : 

D\D,N -K,,I3')^Q . (135) 

Note that this equ atio n is fully general and does not assume that the 3- metric 7 is conformally flat. The vanishing 
of the divergence (|135| ) enables one to use the Green-Ostrogradski formula to get an identity involving only surface 
integrals: 

2 

/ {D^N - K^jl3')dS' = - ^ / {D^N ~ K^jl3^)dS' , (136) 

"'00 0=1 "''5" 

where by convention dS** is oriented towards the "interior" of the throats. 
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From Eq. (13C), the flux integral of D^N on the left hand side is equal to AtiM. Using Eqs. {V2t) and (108), the flux 
integral of is equ al to 87rilJ. The second term on the right hand side vanishes because /3 = on the throats 

[rigidity condition, Eq. (Ill)], so that one is left with 



M -2nj - 



1 

47r 



D,NdS' 



1 

47r 



D.NdS' 



(137) 
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This formula generalizes to the binary black hole case the classical formula that Smarr derived for a single rotating 
black hole (the surface integral on the right hand side being then the black hole surface gravity multiplied by the 
horizon area). 



IV. ASYMPTOTIC BEHAVIOR OF THE FIELDS 



The asy mpto tic b ehavior (near spatial infinity) of the conformal factor ^! and the lapse function N are given by 
Eqs. (129) and ( 130 ). The aim of this section is to get the asymptotic behavior of the shift vector (3 and the extrinsic 
curvature tensor K (or equivalently Aij). In doin g so, we will gain some insight about the assumption of asymptotic 
flatness and the leftover Einstein equations (102). 

To simplify the analysis, we restrict it to a system of identical black holes. We introduce a Cartesian coordinate 
system (a;,y, z) such that x is the direction along the two hole centers (i.e. the centers of the spheres 5*1 and 6*2), 
x = at the middle between the two, and z is the direction perpendicular to the orbital plane. Moreover, let us 
introduce coordinate systems centered on each hole, according to 




and 




(138) 



where d denotes the coordinate distance between the centers of the two spheres Si and iS'2. These two coordinate 
systems are represented in Fig. |4| 




FIG. 4. Cartesian coordinate systems used for the computation of the asymptotic behavior of the shift vector. 



A. Asymptotic behavior of the shift vector 



Let us split Eq. (|100[) for the shift vector in two parts, assuming that its right-hand side can be split in a part S\ 



centered on hole 1 and another part centered ^2 on hole 2. We therefore write 



13' ^9. rd 



2 i 



(139) 
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where is the rotational Killing vector d/dip of the flat metric f already introduced in Sec. Ill C, and (3\ and /?2 are 
the asymptotically vanishing solutions of 



1 



A/3: + 7:D'D,Pi = SI , 



1, 2 . 



Let us solve Eq. ( |140| ) for a = 1 by means of the following decomposition [|65| 

7, 



PI = 



I dxi dwl , 
8 I dx\ dx\ 1 



(140) 



(141) 



where x\ denotes the Cartesian coordinate system (xi,?/i,zi) centered on hole 1, and components with respect to 
that coordinate system are understood in Eq. (141). Wl and xi are solutions of the Poisson equations 



4 



(142) 
(143) 



Provided that the source SI decays at least as r-^ as ri 
is given by the harmonic function 



oo 



o,-,o 

n 



D , the leading behavior of the solution to Eq. ( |142| ) 

+ O(rr') , (144) 



where a is a constant. Note that we have neglected the monopolar part of S^^ and S^^ with respect to that of Sf^ . 
This amounts to considering that S{ corresponds mainly to a motion along the yi direc tion , in accordance with the 
orbital motion. To understand this, let us note that taking the Laplacian of expression ( |l44|) results in the following 
form for S*}: 



(5r,5f ,5^) = (0,-47ra<5(a;i,yi,zi),0) , 
where 6 denotes the Dirac distribution. The Newtonian limit for SI is 

SI = Wirpv' , 



(145) 
(146) 



where p and denotes the matter density and veloc ity r espectively. This Newtonian limit holds because in presence 
of matter the right-hand side of the shift equation (100) should contain the term IGtt times the momentum density 
of matter (see e.g. Eq. (52) of Ref. |^^). For two point mass particles of individual mass m in circular orbit with 
angular velocity O, this results in 



{S^\Sf\S'i') = iO,-87TmdnS{xi,yi,zi),0) 



Identification with (14^) leads to the Newtonian value of the coefficient a: 

ONcwt — 2mdfl . 



(147) 



(148) 



Regarding the Poiss on eq uation (143) for xi, we notice that its source has a vanishing effective mass, at least if 
its leading order is as (|l45| ); consequently, the soluti on x i has no monopolar term in rj~^ and decays as rf^. This 
means that its gradient - which enters in expression (141) for the shift vector - decays as r^^ . Now, in this section, 
we are interested in the behavior of the shift vector up to the order only. Therefore, we discard the solution for 
Xi, writing 



Inserting (144) and (|149|) into (141) yields 



(/3r,/3f ,/5J^) 



/ axiyi 


a 


7^ 


yl] 


aziyA 


\ 8r? 


'8^ 




' 8r? ; 



(149) 



(150) 



ri = \/xY+yi+Zi is the same radial coordinate as that introduced in Sec. II A 1 
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Let us remark that thi s so lution is nothing but one of the three (harmonic) eigenvectors of the operator D ■ L = 
A + 1/2>D{D-) [cf. Eq. (105)] which decay as . This can be seen by comparing with the hst of these harmonic 
vectors provided by 6 Murchadha |^^: the solution (150) above is the item (19.2) of 6 Murchadha's list. Moreover, 



this author has shown that this harmonic vector is g ener ated from linear momentum in the yi direction, in full 
accordance with the analysis performed above [cf. Eq. (145)]. 

At this stage, our solution (150) describes only the linear momentum of hole 1. Since we are considering corotating 
black holes, they must have individual angular momentum (spin), in addition to their linear momentum, although 
neither the notion of individual spin nor individual linear momentum can be defi ned r igorously for a binary system 
in general relativity (only the total angular momentum can be defined, as in Sec. [II C ). To take the rotation of the 
black holes into account, let us add a pure spin part to (3\, of the type 



( B""^ B""^ 



2s 



yi 



-2s^,0 



(151) 



where the constant s is some p aram eter which measures the amount of spin, the latter being supposed to be aligned 
along the zi axis. Note that (151) is a harmonic vector of the operator D ■ L which decays as r^^. It is nothing 
but the asymptotic part of the axisymmetri c sh ift v ector generated by a single rotating object [see e.g. Eq. (4.13) 
of Ref. [|6^, where = —P^pi^]- Adding (15C) to (151), we get the following final expression for the shift vector 
"mostly generated" by hole 1: 



/3f = 



8rf 

a 

aziyi 
8r3 



7- 



2,g + 0(^-3) 



vl 



+ 0(rr^) 



(152) 

(153) 
(154) 



By symmetry, we get exactly the same expression for the components of the shift vector with respect to the 
coordinates (x2, y2, Z2)- Let us now express the components of both fl^ and with respect to the coordinates (x, y, z) 
centered on the system. Taking into account the orientations of (a;i,yi,zi) and (a;2, 2/2,^2) with respect to {x,y,z) 
(see Fig. ^), we obtain: 



PI 

PI 
PI 



a{x + d/2)y 



8rf 



2s4 + 0(r-3) 



(155) 

(156) 
(157) 



and 



P2 
PI 
PI 



a{d/2 — x)y 



8r3 



+ 2s4 + 0(r-3) 



a 

8^ 
azy 

'8^ 



^'^-2,s:^ + 0(.-3) 



+ 0(r-') 



(158) 

(159) 
(160) 



By adding together these two expressions [cf. Eq. ( |139| )] and performing an expansion to the order r ^, we get the 
following asymptotic form of the total shift vector (3: 



py^nx-'^^^^ilx^ + lQy^ 
8 r° 

r = -^^ + 0(r-3). 



7z^)-As^ + 0{r-^) 



(161) 
(162) 
(163) 
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Note that, apart from the part, the total shift decays as r contrary to (3i and /Sj, which decay as r ^ . From the 
above expression, /3 can be hnearly decomposed into three parts: 

/3 = /3kin + /3angu + /^quad > (164) 

with 

Ptin = -^y (165) 

fil^ = nx (166) 

/3kin = , (167) 

/3aV--(T+4^)^ (169) 

/3a„gu = (170) 

and 



f}X _ 

Hquad 


iad V ( 

H 


a;2 

" ,.2 


Ry - 

r^quad 


3ad X ( 
- 8 r^\^' 


-4 


— 

Hquad 


3ad xyz 
8 r5 





(171) 
(172) 
(173) 

/3]jijj is a pure kinematical term, which reflects that we use co-rotating coordinates. /Q^j^gu is, as we will see belo w, th e 
part of the shift which carries the total angular-momentum of the system. As for ^^-^^ introduced above [Eq. 151)], 
it has the familiar shape of a pure spin axisymmetric shift vector. 

/3quad is oil's of thc nine harmonic vectors of the operator D ■ L which decay as r^"^ . It has the number (22.7) in O 
Murchadha's list By the way, /3angu the number (22.1) in the same hst. As shown by 6 Murchadha |Q [cf. 
his Eq. (29)], /3quad arises from the fact that the Qxy component of the quadrupole moment Qij of the system is time 
varying with respect to some asymptotic inertial frame ^. This is the only such component. Indeed, if we consider a 
Newtonian system of two identical point mass particles on a circular orbit of diameter d in the x — y plane, the time 
derivative of its quadrupole moment with respect to thc inertial frame is given by 

^ -m—nsm{2nt) (174) 
2 

d^ 
2 

^Qzz^O. (177) 



Qxy ^ m—ncos{2nt) (175) 
Qyy = m^nsm{2nt) (176) 



It is clear on this expression that at time i = 0, when the axes of the rotating frame coincide with those of the 
inertial frame (our assumption in this discussion), the only non- vanishing component is Q^y = md^il/2. Combining 



Eqs. (|148| ) and (|l75| ), we see that the coefficient in front of the three components (|l7l| )-( |l73[) of /Sq^^^^j is — 3/2Qx 



This justifies the label quad (for quadrupole moment) given to that part of the shift vector. 



^In Post-Newtonian theory, it is also well known that some part of the gravitomagnetic potential - the shift vector in our 
language - can be generated by the first time derivative of the mass quadrupole moment, see e.g. Sec. VLB of Ref. jHsj or 
Eq. (4.2) of Ref. m. 
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B. Asymptotic behavior of the extrinsic curvature 



Th e asymptotic behavior of the extrinsic curvature tensor is deduced from that of the shift vector via Eqs. (|9^) and 
( |104| ), which aUows us to write (taking into account that both and N are equal to 1 at spatial infinity) 



when r —>■ oo 



(178) 



Since (L/?kin)'-' — for /J]^;^ is a Killing vector of the flat metric f, the decomposition (164) of the shift vector leads 
to the following decomposition of the extrinsic curvature tensor: 



K = K 



angu 



K 



quad 



(179) 



where -?^angu i^PanguY^ ^nd i^quad i^Pquad)^^ • Inserting the formulas (168)-(I70) in the exphcit form (105) of 
the operator L results in the following components of Kangu^ 



angu 



angu 



angu 



Lryy 

angu 



angu 
angu 



3 ( — + 2s 



9 9 

X — y 



= -3(- + 2s|^ 



ad 



= 3 — + 4s ^ 



xy 



„ , ad „ , _.- 
= 3 — + 2s — 



= 



(180) 

(181) 

(182) 

(183) 

(184) 
(185) 



A similar computation for K^^^^ yields 



ly-xx 
quad 



quad 



T^XZ 

quad 



T^VV 
quad 



quad 



K. 



quad 



^adxyfx^^^ 



= 3 



9 9 
5 Z 



8r5 

y 

8 rS 



3^^ ("5^ + 1 



8 r5 



= 3—^ 5^ + 1 



= 3 



8 r 



adxy f z 



5^-3 



(186) 
(187) 
(188) 
(189) 
(190) 
(191) 



Note that both K^i^ and ifqu^d decay as r ^. 



If we plug the formulas (18(;)-(191) into the surface integral (125) which gives the angular momentum, we get, after 
some straightforward calculations 



'^(Kquad) — 



(192) 



This means th at all the angular momentum of the system is carried by K'^i^^^. Indeed, inserting the formulas ( |180|) - 
(|T8|) into Eq. (|l2|) gives 



T/., N Old 

J = J(Kangu) = — + 2S 



(193) 
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For a non-relativistic point mass system, a is given by Eq. (148) so that we get 



^(Ka„gu)Ncwt = m—Q. + 2s . (194) 

The first term on the right-hand side is nothing but the orbital angular momentum of the system and the second 
terms is the sum of the spins s of the two particles. Hence J(Kangu)Newt is equal to the total angular momentum of 
the system. 

C. Helical symmetry and asymptotic flatness 

Let us consider the five Einstein equations that we have not taken into account for the solution of the problem, i.e. 



Eqs. (102). Thanks to the asymptotic behavior ( |129| ), ( |130| ) and (131) of N and all the terms involving products 
of gradients of N or ^, as well as the ones involving second derivatives of ^^iV, decay at least as r~^. The quadratic 



term fkiA^^A'^ decays as r ^ for K"^^ decays as r ^ , as seen above. The only remaining term in Eq. (102) is the Lie 
derivative of A^^ along /3. Asymptotically, one has 

£fjA'^ = £^K'' + 0{r-^) . (195) 

It can be seen easily that only the Lie derivative along /Sj^jj, matters: 

£fjK'^ = £fj^ K'^ + 0{r~'^) . (196) 



Let us introduce the splitting (17£) of K^^ into this expression. After some computations, we find that 

0(r-^) , (197) 



whereas 



£(3K^..i = ^(42^ - 5x^r - Sx^z^ + - v'z' - 2z') + ©(r"^) (198) 

£pKL. = - V') + 0{r-') (199) 

-^/S^quad = ^2;z(7x2 - + 2z2) + 0(^-4) (200) 

= ^(-^' + + - V + 32;2z^ + 2z4) + 0(,-4) (201) 

£(3KL^ - ^-^y<^-' - - 2.^) + 0{r-^) (202) 

"^/3^quad - -^(-' - y')(3x^ + 3,/ - 2z^) + 0{r~^) , (203) 

which means that £ j^K^^^^ decays only as . We face here the incompatibility of helical symmetry and asymptotic 
flatness for systems that have a time-varying quadrupole moment (recall that /3quad hence Kquad is due to 
Qxy 7^ 0). Indeed £pK^^^^ is the only term in the Einstein equations (102) which decays as slower as r~^. It 



therefore cannot be compensated by another term. This means that the five Einstein equations ( |102| ) are violated. 
Note that this problem does not arise from the assumption of conformal flatness of the 3-metric 7. Relaxing this 
condition would have resulted in asymptotic behaviors of f3 and K which would have been the same as that obtained 
here. Note also that for a system such as an isolated rotating axisymmetric star (or m ore g enerally for any stationary 
system), /Sq^ad ~ and Kquad = 0, so that the problem of asymptotic flatness in Eq. (102) does not arise. 



V. CONCLUSIONS 

We have presented an approach to the problem of binary black holes in circular orbit which is similar to that 
previously used in the literature to treat binary neutron stars [p2|-p9| , namely an approach based on the existence 



24 



of a helical Killing vector field along with the simplifying assumption of a conformally flat 3-metric. The differences 
between the two approaches lie in the boundary conditions on the throats in the black hole case. We have chosen a 
spacetime manifold with spatial sections of the Misner-Lindquist type, i.e. composed of two isometric asymptotically 
flat sheets. Moreover, we have chosen the black holes to be corotating, which has the simple geometrical interpretation 
of the throats being Killing horizons. 

By enforcing the isometry conditions on the shift vector, as well as the equation relating the trace of the extrinsic 
curvature to the divergence of the shift, possibly at the price of slightly modifying the momentum constraint, all the 
quantities which enter in the equations remain regular. Notably the extrinsic curvature tensor remains finite on the 
throats, although the lapse vanishes there. 

We have proposed to compute the orbital angular velocity of the system by requiring that the conformal factor 
^' and the lapse function N have the same monopolar 1/r term in their asymptotic expansions. This requirement 
reduces to the classical virial theorem at the Newtonian limit. 

Contrary to previous numerical approaches - the conformal imaging one p^[2^ and the puncture one |2^j2^ - our 
method amounts to solving five, and not four (the four constraints), of the Einstein equations. This reflects the fact 
that we have re-introduced the time dimension in the problem. 

The formulation presented here has been implemented by means of a numerical code based on a multi-domain 
spectral method and we present the first results in the companion paper [ pO[ . These results can be used as initial 
data for computing the black hole coalescence within the 3-1-1 formalism [p|-|ll|. 

Let us stress that the work presented in this article constitutes a first attempt to tackle the problem of binary 
black hole in circular orbits. In order to fully specify the problem and search for a unique solution, we had to make a 
number of concrete choices which have some degree of arbitrariness, such as the two-sheeted topology, the isometry 
across the throats and the resulting boundary conditions, or the rigid rotation of the black holes. These hypotheses 
could be changed to different ones, as for instance considering irrotational black holes instead of corotating ones. This 
shall be investigated in future works. 
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